home *** CD-ROM | disk | FTP | other *** search
/ Amiga Plus 1995 #5 & #6 / Amiga Plus CD - 1995 - No. 5 and 6.iso / pd / daten / astrolog / src / placalc2.c < prev    next >
Encoding:
C/C++ Source or Header  |  1995-08-11  |  35.8 KB  |  956 lines
  1. /*                                                               -*- C -*-
  2. ** Astrolog (Version 4.40) File: placalc2.c
  3. **
  4. ** IMPORTANT NOTICE: The graphics database and chart display routines
  5. ** used in this program are Copyright (C) 1991-1995 by Walter D. Pullen
  6. ** (astara@u.washington.edu). Permission is granted to freely use and
  7. ** distribute these routines provided one doesn't sell, restrict, or
  8. ** profit from them in any way. Modification is allowed provided these
  9. ** notices remain with any altered or edited versions of the program.
  10. **
  11. ** The main planetary calculation routines used in this program have
  12. ** been Copyrighted and the core of this program is basically a
  13. ** conversion to C of the routines created by James Neely as listed in
  14. ** Michael Erlewine's 'Manual of Computer Programming for Astrologers',
  15. ** available from Matrix Software. The copyright gives us permission to
  16. ** use the routines for personal use but not to sell them or profit from
  17. ** them in any way.
  18. **
  19. ** The PostScript code within the core graphics routines are programmed
  20. ** and Copyright (C) 1992-1993 by Brian D. Willoughby
  21. ** (brianw@sounds.wa.com). Conditions are identical to those above.
  22. **
  23. ** The extended accurate ephemeris databases and formulas are from the
  24. ** calculation routines in the program "Placalc" and are programmed and
  25. ** Copyright (C) 1989,1991,1993 by Astrodienst AG and Alois Treindl
  26. ** (alois@azur.ch). The use of that source code is subject to
  27. ** regulations made by Astrodienst Zurich, and the code is not in the
  28. ** public domain. This copyright notice must not be changed or removed
  29. ** by any user of this program.
  30. **
  31. ** Initial programming 8/28,30, 9/10,13,16,20,23, 10/3,6,7, 11/7,10,21/1991.
  32. ** X Window graphics initially programmed 10/23-29/1991.
  33. ** PostScript graphics initially programmed 11/29-30/1992.
  34. ** Last code change made 1/29/1995.
  35. */
  36.  
  37. #include "placalc.h"
  38.  
  39.  
  40. #ifdef PLACALC
  41. /*
  42. ** ---------------------------------------------------------------
  43. ** | Copyright Astrodienst AG and Alois Treindl, 1989,1991,1993  |
  44. ** | The use of this source code is subject to regulations made  |
  45. ** | by Astrodienst Zurich. The code is NOT in the public domain.|
  46. ** |                                                             |
  47. ** | This copyright notice must not be changed or removed        |
  48. ** | by any user of this program.                                |
  49. ** ---------------------------------------------------------------
  50. **
  51. ** Important changes:
  52. ** 11-jun-93 revision 1.12: fixed error which affected Mercury between -2100
  53. ** and -3100 (it jumped wildly).
  54. */
  55.  
  56. /***********************************************************
  57. ** $Header$
  58. **
  59. ** definition module for planetary elements
  60. ** and disturbation coefficients
  61. ** version HP-UX C  for new version with stored outer planets
  62. ** 31-jul-88
  63. ** by Alois Treindl 
  64. **
  65. ** ---------------------------------------------------------------
  66. ** | Copyright Astrodienst Zurich AG and Alois Treindl, 1989.    |
  67. ** | The use of this source code is subject to regulations made |
  68. ** | by Astrodienst Zurich. The code is NOT in the public domain.|
  69. ** |                                                             |
  70. ** | This copyright notice must not be changed or removed        |
  71. ** | by any user of this program.                                |
  72. ** ---------------------------------------------------------------
  73. **
  74. ***********************************************************/
  75.  
  76. /*
  77. ** In the elements degrees were kept as the units for the constants. This
  78. ** requires conversion to radians, when the actual calculations are performed.
  79. ** This approach is not the most efficient, but safer for development.
  80. ** Constant conversion could be done by writing all degree constants with
  81. ** value * DEGTORAD
  82. */
  83.  
  84. #define TIDAL_26  TRUE  /* decide wheter to use new or old lunar tidal
  85. term; a consistent system of delta t must be 
  86. used */
  87. #define MOON_TEST_CORR FALSE  /* to include more lunar terms in longitude */
  88.  
  89. REAL8 ekld[4] = {23.452294, -46.845, -.0059, 0.00181};
  90. /* ecliptic with epoch1900, Ekd(0..3) in basic */
  91.  
  92. struct eledata pd[MARS + 1] = {
  93.  
  94.   {  /* earth */
  95.  
  96.         1.00000023,  365.25636042,    EPOCH1900,
  97.        99.696678,      0.9856473354,  1.089,         0,
  98.       101.220833,   6189.03,          1.63,          0.012,
  99.         0.01675104,   -0.00004180,   -0.000000126,
  100.         0,             0,             0,    0,
  101.         0,             0,             0
  102.   },
  103.  
  104.   /*
  105.   ** note 29 June 88 by Alois: G.M.Clemence, Astronomical Journal
  106.   ** vol.53,p. 178 (1948) gives a correction to the perihel motion
  107.   ** of -4.78" T, giving 6184.25 for the linear Term above. We have
  108.   ** not yet applied this correction. It has been used in APAE 22,4
  109.   ** on the motion of mars and does make an official impression.
  110.   */
  111.  
  112.   {  /* moon */ 
  113.  
  114.     0.0025955307, 27.321661,  EPOCH1900,
  115.  
  116. #if ! TIDAL_26
  117. /*
  118. ** values from Improved Lunar Ephemeris, corresponding to tidal
  119. ** term -22.44"/cy and  consistent with delta t ~ 29.949 T*T
  120. */
  121.  
  122.     270.4341638,  13.176396526808121, -4.08,  0.0068,
  123. #endif
  124.  
  125.  
  126. #if TIDAL_26
  127. /*
  128. ** new values from Morrison 1979, with tidal term -26"/cy as
  129. ** stated in A.E. 1986 onwards, consistent with delta t ~ 44.3 T*T
  130. ** correction: -1.54" + 2.33" T - 1.78" T*T
  131. */
  132.  
  133.     270.4337361,  13.176396544528099, -5.86,  0.0068,
  134. # endif
  135.  
  136.  
  137.     334.329556, 14648522.52,  -37.17,   -0.045,
  138.     0.054900489,  0,    0,
  139.     259.183275, -6962911.23,  7.48 ,    0.008,
  140.     5.145388889,  0,    0
  141.   },
  142.  
  143.   { /* mercury */
  144.  
  145.     0.3870986, 87.969252,  EPOCH1900,
  146.     178.179078, 4.0923770233, 1.084,    0,
  147.     75.89969722,  5599.76,  1.061,    0,
  148.     0.20561421, .00002046,   -.000000030,
  149.     47.145944444, 4266.75,  .626,   0,
  150.     7.0028805555, 6.699,    -.066
  151.   },
  152.  
  153.  
  154.   { /* venus */
  155.  
  156.     0.72333162,  224.700726, EPOCH1900,
  157.     342.767053, 1.6021687039, 1.1148,   0,
  158.     130.16383333, 5068.93,  -3.515,   0,
  159.     0.00682069, -.00004774, .000000091,
  160.     75.7796472223,3239.46,  1.476,    0,
  161.     3.3936305555, 3.621,    .0035
  162.   },
  163.  
  164.  
  165.   { /* mars */ 
  166.  
  167.     1.5236914620, 686.9296097,  EPOCH1900,
  168.  
  169.   /* These are the corrected elements by Ross */
  170.  
  171.     293.74762778, .524071163814,  1.1184,   0,
  172.     334.21820278, 6626.73,  .4675,    -0.0043,
  173.     0.09331290, .000092064, -.000000077,
  174.     48.786441667, 2775.57,  -.005,    -0.0192,
  175.     1.85033333, -2.430,   .0454
  176.   }
  177. };
  178.  
  179. /*
  180. ** mimimum and maximum distances computed over 1000 years with plamimax,
  181. ** required for relative distances rgeo, where the distance is given
  182. ** as 100 when a planet is closest and as 0 when farthest from earth.
  183. */
  184.  
  185. REAL8 rmima[CALC_N][2] = {
  186.   {   0.98296342,    1.01704665 },
  187.   {   0.00238267,    0.00271861 },
  188.   {   0.54900496,    1.45169607 },
  189.   {   0.26411287,    1.73597885 },
  190.   {   0.37289847,    2.67626927 },
  191.   {   3.94877993,    6.45627627 },
  192.   {   7.99362824,   11.09276636 },
  193.   {  17.28622633,   21.10714104 },
  194.   {  28.81374786,   31.33507284 },
  195.   {  28.67716748,   50.29208774 },
  196.   {   0.00,          0.00259553 }, /* nodes don't get a real value */
  197.   {   0.00,          0.00259553 },
  198.   {   7.36277475,   19.86585062 }
  199. };
  200.  
  201.  
  202. struct sdat _sd [SDNUM] = {
  203.   { 114.50,       585.17493     },
  204.   { 109.856,      191.39977     },
  205.   { 148.031,       30.34583     },
  206.   { 284.716,       12.21794     },
  207.   { 114.508,      585.17656     },
  208.   {  -0.56,       359.99213     },
  209.   { 148.03,        30.34743     },
  210.   { 284.72,        12.2196      },
  211.   { 248.07,      1494.726615    },
  212.   { 359.44,       359.993595    },
  213.   { 109.86,       191.402867    },
  214.   { 148.02,        30.348930    },
  215.   { 114.503,      585.173715    },
  216.   { 359.444,      359.989285    },
  217.   { 148.021,       30.344620    },
  218.   { 284.716,       12.21669     },
  219.   { 148.0315,      30.34906264  },
  220.   { 284.7158,      12.22117085  },
  221.   { 220.1695,       4.284931111 },
  222.   { 291.8024,       2.184704167 }
  223. };
  224.  
  225. REAL8 sa[SDNUM];
  226.  
  227. /*
  228. ** delta long = lampl * COS (lphase - arg) in seconds of arc
  229. ** delta rad  = rampl * COS (rphase - arg) in ninth place of log
  230. ** arg = j * sa (k) + i * ma (this planet)
  231. ** ma = mean anomaly
  232. ** sa = mean anomaly of disturbing planet, where this
  233. ** is taken from the aproximate value in sa[]
  234. ** For the COS (phase - arg) it is good enough to compute
  235. ** with 32 bit reals, because ampl and phase have only 
  236. ** four to five significant digits.
  237. ** While saving constant space, it is costing execution time due
  238. ** to float/double conversions.
  239. **
  240. ** In basic, all correction terms for sun, mercury, venus and mars
  241. ** were contained in one array K(0..142,0..6); Nk(N,0) contained
  242. ** the index of the first term of planet N and Nk(N,1) the number
  243. ** of terms for this planet. Here, we use a  0 in the first column
  244. ** kor.j to indicate the end of the table for a planet.
  245. ** K(*) was a basic INTEGER array, therefore the amplitudes and phases
  246. ** had to be expressed as
  247. ** K(i,2) = ampl. of longitude in 0.001 seconds of arc
  248. ** K(i,3) = phase of longitude in 0.01 degrees
  249. ** K(i,4) = ampl. of radius in 9th place of log
  250. ** K(i,5) = phase of radius in 0.01 degrees.
  251. ** Here we have converted the amplitude of long. to seconds of arc
  252. ** and the phases to degrees.
  253. */
  254.  
  255. struct kor ARR earthkor[] = {  /* 11-jul-88 all terms to 0.020" longitude */
  256.  
  257.   /* j   i    lampl   lphase  rampl   rphase    k */
  258.  
  259.   { -1,  1,   0.013,  243,      28,   335,      8 },  /* mercury */
  260.   { -1,  3,   0.015,  357,      18,   267,      8 },
  261.   { -1,  4,   0.023,  326,       5,   239,      8 },
  262.   { -1,  0,   0.075,  296.6,    94,   205.0,    0 },  /* venus */
  263.   { -1,  1,   4.838,  299.10, 2359,   209.08,   0 },
  264.   { -1,  2,   0.074,  207.9,    69,   348.5,    0 },
  265.   { -1,  3,   0.009,  249,      16,   330,      0 },
  266.   { -2,  1,   0.116,  148.90,  160,    58.40,   0 },
  267.   { -2,  2,   5.526,  148.31, 6842,    58.32,   0 },
  268.   { -2,  3,   2.497,  315.94,  869,   226.70,   0 },
  269.   { -2,  4,   0.044,  311.4,    52,    38.8,    0 },
  270.   { -3,  2,   0.013,  176,      21,    90,      0 },
  271.   { -3,  3,   0.666,  177.71, 1045,    87.57,   0 },
  272.   { -3,  4,   1.559,  -14.75, 1497,   255.25,   0 },
  273.   { -3,  5,   1.024,  318.15,  194,    49.50,   0 },
  274.   { -3,  6,   0.017,  315,      19,    43,      0 },
  275.   { -4,  4,   0.210,  206.20,  376,   116.28,   0 },
  276.   { -4,  5,   0.144,  195.40,  196,   105.20,   0 },
  277.   { -4,  6,   0.152,  -16.20,   94,   254.80,   0 },
  278.   { -5,  5,   0.084,  235.6,   163,   145.4,    0 },
  279.   { -5,  6,   0.037,  221.8,    59,   132.2,    0 },
  280.   { -5,  7,   0.123,  195.30,  141,   105.40,   0 },
  281.   { -5,  8,   0.154,    -.40,   26,   270.00,   0 },
  282.   { -6,  6,   0.038,  264.1,    80,   174.3,    0 },
  283.   { -6,  7,   0.014,  253,      25,   164,      0 },
  284.   { -6,  8,   0.01,   230,      14,   135,      0 },
  285.   { -6,  9,   0.014,   12,      12,   284,      0 },
  286.   { -7,  7,   0.020,  294,      42,   203.5,    0 },
  287.   { -7,  8,   0.006,  279,      12,   194,      0 },
  288.   { -8,  8,   0.011,  322,      24,   234,      0 },
  289.   { -8, 12,   0.042,  259.2,    44,   169.7,    0 },
  290.   { -8, 14,   0.032,   48.8,    33,   138.7,    0 },
  291.   { -9,  9,   0.006,  351,      13,   261,      0 },
  292.   {  1, -1,   0.273,  217.70,  150,   127.70,   1 }, /* mars */
  293.   {  1,  0,   0.048,  260.3,    28,   347,      1 },
  294.   {  2, -3,   0.041,  346,      52,   255.4,    1 },
  295.   {  2, -2,   2.043,  343.89, 2057,   253.83,   1 },
  296.   {  2, -1,   1.770,  200.40,  151,   295.00,   1 },
  297.   {  2,  0,   0.028,  148,      31,   234.3,    1 },
  298.   {  3, -3,    .129,  294.20,  168,   203.50,   1 },
  299.   {  3, -2,    .425,  -21.12,  215,   249.00,   1 },
  300.   {  4, -4,   0.034,   71,      49,   339.7,    1 },
  301.   {  4, -3,    .500,  105.18,  478,    15.17,   1 },
  302.   {  4, -2,    .585,  -25.94,  105,    65.90,   1 },
  303.   {  5, -4,   0.085,   54.6,   107,   324.6,    1 },
  304.   {  5, -3,    .204,  100.80,   89,    11.00,   1 },
  305.   {  6, -5,   0.020,  186,      30,    95.7,    1 },
  306.   {  6, -4,    .154,  227.40,  139,   137.30,   1 },
  307.   {  6, -3,    .101,   96.30,   27,   188.00,   1 },
  308.   {  7, -5,   0.049,  176.5,    60,    86.2,    1 },
  309.   {  7, -4,    .106,  222.70,   38,   132.90,   1 },
  310.   {  8, -5,   0.052,  348.9,    45,   259.7,    1 },
  311.   {  8, -4,   0.021,  215.2,     8,   310,      1 },
  312.   {  8, -6,   0.010,  307,      15,   217,      1 },
  313.   {  9, -6,   0.028,  298,      34,   208.1,    1 },
  314.   {  9, -5,   0.062,  346,      17,   257,      1 },
  315.   { 10, -6,   0.019,  111,      15,    23,      1 },
  316.   { 11, -7,   0.017,   59,      20,   330,      1 },
  317.   { 11, -6,   0.044,  105.9,     9,    21,      1 },
  318.   { 13, -8,   0.013,  184,      15,    94,      1 },
  319.   { 13, -7,   0.045,  227.8,     5,   143,      1 },
  320.   { 15, -9,   0.021,  309,      22,   220,      1 },
  321.   { 17, -9,   0.026,  113,       0,     0,      1 },
  322.   {  1, -2,    .163,  198.60,  208,   112.00,   2 }, /* jupiter */
  323.   {  1, -1,   7.208,  179.53, 7067,    89.55,   2 },
  324.   {  1,  0,   2.600,  263.22,  244,   -21.40,   2 },
  325.   {  1,  1,   0.073,  276.3,    80,     6.5,    2 },
  326.   {  2, -3,   0.069,   80.8,   103,   350.5,    2 },
  327.   {  2, -2,   2.731,   87.15, 4026,    -2.89,   2 },
  328.   {  2, -1,   1.610,  109.49, 1459,    19.47,   2 },
  329.   {  2,  0,   0.073,  252.6,     8,   263,      2 },
  330.   {  3, -3,    .164, 170.50,   281,    81.20,   2 },
  331.   {  3, -2,    .556,  82.65,   803,    -7.44,   2 },
  332.   {  3, -1,    .210,  98.50,   174,     8.60,   2 },
  333.   {  4, -4,   0.016, 259,       29,   170,      2 },
  334.   {  4, -3,   0.044, 168.2,     74,    79.9,    2 },
  335.   {  4, -2,   0.080,  77.7,    113,   347.7,    2 },
  336.   {  4, -1,   0.023,  93,       17,     3,      2 },
  337.   {  5, -2,   0.009,  71,       14,   343,      2 },
  338.   {  1, -2,   0.011, 105,       15,    11,      3 }, /* saturn */
  339.   {  1, -1,    .419, 100.58,   429,    10.60,   3 },
  340.   {  1,  0,    .320, 269.46,     8,    -7.00,   3 },
  341.   {  2, -2,    .108, 290.60,   162,   200.60,   3 },
  342.   {  2, -1,    .112, 293.60,   112,   203.10,   3 },
  343.   {  3, -2,   0.021, 289,       32,   200.1,    3 },
  344.   {  3, -1,   0.017, 291,       17,   201,      3 },
  345.  
  346.   {ENDMARK}
  347. };
  348.  
  349. struct kor ARR mercurykor[] = {
  350.  
  351.   /* j   i    lampl   lphase  rampl   rphase    k */
  352.  
  353.   {  1,  -1,   .711,  35.47,   491,   305.28,   4 },
  354.   {  2,  -3,   .552, 161.15,   712,    71.12,   4 },
  355.   {  2,  -2,  2.100, 161.15,  2370,    71.19,   4 },
  356.   {  2,  -1,  3.724, 160.69,   899,    70.49,   4 },
  357.   {  2,   0,   .729, 159.76,   763,   250.00,   4 },
  358.   {  3,  -3,   .431, 105.37,   541,    15.53,   4 },
  359.   {  3,  -2,  1.329, 104.78,  1157,    14.84,   4 },
  360.   {  3,  -1,   .539, 278.95,    14,   282.00,   4 },
  361.   {  4,  -2,   .484, 226.40,   234,   136.02,   4 },
  362.   {  5,  -4,   .685, -10.43,   849,   259.51,   4 },
  363.   {  5,  -3,  2.810, -10.14,  2954,   259.92,   4 },
  364.   {  5,  -2,  7.356, -12.22,   282,   255.43,   4 },
  365.   {  5,  -1,  1.471, -12.30,  1550,    77.75,   4 },
  366.   {  5,   0,   .375, -12.29,   472,    77.70,   4 },
  367.   {  2,  -1,   .443, 218.48,   256,   128.36,   5 },
  368.   {  4,  -2,   .374, 151.81,   397,    61.63,   5 },
  369.   {  4,  -1,   .808, 145.93,    13,    35.00,   5 },
  370.   {  1,  -1,   .697, 181.07,   708,    91.38,   6 },
  371.   {  1,   0,   .574, 236.72,    75,   265.40,   6 },
  372.   {  2,  -2,   .938,  36.98,  1185,   306.97,   6 },
  373.   {  2,  -1,  3.275,  37.00,  3268,   306.99,   6 },
  374.   {  2,   0,   .499,  31.91,   371,   126.90,   6 },
  375.   {  3,  -1,   .353,  25.84,   347,   295.76,   6 },
  376.   {  2,  -1,   .380, 239.87,     0,     0,      7 },
  377.  
  378.   {ENDMARK}
  379. };
  380.  
  381. struct kor ARR venuskor[] = {
  382.  
  383.   /* j   i    lampl   lphase  rampl   rphase    k */
  384.  
  385.   { -1,   2,   .264, -19.20,   175,   251.10,   8 },
  386.   { -2,   5,   .361, 167.68,    55,    77.20,   8 },
  387.   {  1,  -1,  4.889, 119.11,  2246,    29.11,   9 },
  388.   {  2,  -2, 11.261, 148.23,  9772,    58.21,   9 },
  389.   {  3,  -3,  7.128,  -2.57,  8271,   267.42,   9 },
  390.   {  3,  -2,  3.446, 135.91,   737,    47.37,   9 },
  391.   {  4,  -4,  1.034,  26.54,  1426,   296.49,   9 },
  392.   {  4,  -3,   .677, 165.32,   445,    75.70,   9 },
  393.   {  5,  -5,   .330,  56.88,   510,   -33.36,   9 },
  394.   {  5,  -4,  1.575, 193.93,  1572,   104.21,   9 },
  395.   {  5,  -3,  1.439, 138.08,   162,   229.90,   9 },
  396.   {  6,  -6,   .143,  84.40,   236,    -5.80,   9 },
  397.   {  6,  -5,   .205,  44.20,   256,   314.20,   9 },
  398.   {  6,  -4,   .176, 164.30,    70,    75.70,   9 },
  399.   {  8,  -5,   .231, 180.00,    25,    75.00,   9 },
  400.   {  3,  -2,   .673, 221.62,   717,   131.60,  10 },
  401.   {  3,  -1,  1.208, 237.57,    29,   149.00,  10 },
  402.   {  1,  -1,  2.966, 208.09,  2991,   118.09,  11 },
  403.   {  1,   0,  1.563, 268.31,    91,    -7.60,  11 },
  404.   {  2,  -2,   .889, 145.16,  1335,    55.17,  11 },
  405.   {  2,  -1,   .480, 171.01,   464,    80.95,  11 },
  406.   {  3,  -2,   .169, 144.20,   250,    54.00,  11 },
  407.  
  408.   {ENDMARK}
  409. };
  410.  
  411. struct kor ARR marskor[] = {
  412.  
  413.   /* j   i    lampl   lphase  rampl   rphase    k */
  414.  
  415.   { -1,   1,   .115,  65.84,   684,   156.14,  12 },
  416.   { -1,   2,   .623, 246.03,   812,   155.77,  12 },
  417.   { -1,   3,  6.368,  57.60,   556,   -32.06,  12 },
  418.   { -1,   4,   .588,  57.24,   616,   147.28,  12 },
  419.   { -2,   5,   .138,  39.18,   157,   309.39,  12 },
  420.   { -2,   6,   .459, 217.58,    82,   128.10,  12 },
  421.   { -1,  -1,   .106,  33.60,   141,   303.45,  13 },
  422.   { -1,   0,   .873,  34.34,  1112,   304.05,  13 },
  423.   { -1,   1,  8.559,  35.10,  6947,   304.45,  13 },
  424.   { -1,   2, 13.966,  20.50,  2875,   113.20,  13 },
  425.   { -1,   3,  1.487,  22.18,  1619,   112.38,  13 },
  426.   { -1,   4,   .175,  22.46,   225,   112.15,  13 },
  427.   { -2,   2,   .150,  18.96,   484,   266.42,  13 },
  428.   { -2,   3,  7.355, 158.64,  6412,    68.62,  13 },
  429.   { -2,   4,  4.905, 154.09,  1985,   244.70,  13 },
  430.   { -2,   5,   .489, 154.39,   543,   244.50,  13 },
  431.   { -3,   3,   .216, 111.06,   389,    21.10,  13 },
  432.   { -3,   4,   .355, 110.64,   587,    19.17,  13 },
  433.   { -3,   5,  2.641, 280.58,  2038,   190.60,  13 },
  434.   { -3,   6,   .970, 276.06,   587,     6.75,  13 },
  435.   { -3,   7,   .100, 276.20,   116,     6.40,  13 },
  436.   { -4,   5,   .152, 232.48,   259,   142.60,  13 },
  437.   { -4,   6,   .264, 230.47,   387,   139.75,  13 },
  438.   { -4,   7,  1.156,  41.64,   749,   312.67,  13 },
  439.   { -4,   8,   .259,  37.92,   205,   128.80,  13 },
  440.   { -5,   8,   .172,  -8.99,   234,   260.70,  13 },
  441.   { -5,   9,   .575, 164.48,   308,    74.60,  13 },
  442.   { -6,  10,   .115, 113.70,   145,    23.53,  13 },
  443.   { -6,  11,   .363, 285.69,   144,   196.00,  13 },
  444.   { -7,  13,   .353,  48.83,    85,   319.10,  13 },
  445.   { -8,  15,  1.553, 170.14,   110,    81.00,  13 },
  446.   { -8,  16,   .148, 170.74,   154,   259.94,  13 },
  447.   { -9,  17,   .193, 293.70,    23,    22.80,  13 },
  448.   {  1,  -3,   .382,  46.48,   521,   316.25,  14 },
  449.   {  1,  -2,  3.144,  46.78,  3894,   316.39,  14 },
  450.   {  1,  -1, 25.384,  48.96, 23116,   318.87,  14 },
  451.   {  1,   0,  3.732, -17.62,  1525,   117.81,  14 },
  452.   {  1,   1,   .474, -34.60,   531,    59.67,  14 },
  453.   {  2,  -4,   .265, 192.88,   396,   103.12,  14 },
  454.   {  2,  -3,  2.108, 192.72,  3042,   102.89,  14 },
  455.   {  2,  -2, 16.035, 191.90, 22144,   101.99,  14 },
  456.   {  2,  -1, 21.869, 188.35, 16624,    98.33,  14 },
  457.   {  2,   0,  1.461, 189.66,  1478,   279.04,  14 },
  458.   {  2,   1,   .167, 191.04,   224,   280.81,  14 },
  459.   {  3,  -4,   .206, 167.11,   338,    76.13,  14 },
  460.   {  3,  -3,  1.309, 168.27,  2141,    76.24,  14 },
  461.   {  3,  -2,  2.607, 228.41,  3437,   139.74,  14 },
  462.   {  3,  -1,  3.174, 207.20,  1915,   115.83,  14 },
  463.   {  3,   0,   .232, 207.78,   240,   298.06,  14 },
  464.   {  4,  -4,   .178, 127.25,   322,    36.16,  14 },
  465.   {  4,  -3,   .241, 200.69,   389,   110.02,  14 },
  466.   {  4,  -2,   .330, 267.57,   413,   179.86,  14 },
  467.   {  4,  -1,   .416, 221.88,   184,   128.17,  14 },
  468.   {  1,  -2,   .155, -38.20,   191,   231.58,  15 },
  469.   {  1,  -1,  1.351, -34.10,  1345,   235.85,  15 },
  470.   {  1,   0,   .884, 288.05,   111,    39.90,  15 },
  471.   {  1,   1,   .132, 284.88,   144,    15.67,  15 },
  472.   {  2,  -2,   .620,  35.15,   869,   305.30,  15 },
  473.   {  2,  -1,  1.768,  32.50,  1661,   302.51,  15 },
  474.   {  2,   0,   .125,  18.73,   103,   119.90,  15 },
  475.   {  3,  -2,   .141,  47.59,   199,   318.06,  15 },
  476.   {  3,  -1,   .281,  40.95,   248,   310.75,  15 },
  477.  
  478.   {ENDMARK}
  479. };
  480.  
  481. struct m45dat m45[NUM_MOON_CORR] = {
  482.  
  483.   /* l    l'   F    D,      Long,       Lat       Par    */
  484.  
  485.    { 0,   0,   0,   4,     13.902,     14.06,    0.2607 },
  486.    { 0,   0,   0,   2,   2369.912,   2373.36,   28.2333 },
  487.    { 1,   0,   0,   4,      1.979,      6.98,    0.0433 },
  488.    { 1,   0,   0,   2,    191.953,    192.72,    3.0861 },
  489.    { 1,   0,   0,   0,  22639.500,  22609.1,   186.5398 },
  490.    { 1,   0,   0,  -2,  -4586.465,  -4578.13,   34.3117 },
  491.    { 1,   0,   0,  -4,    -38.428,    -38.64,    0.6008 },
  492.    { 1,   0,   0,  -6,     -0.393,     -1.43,    0.0086 },
  493.    { 0,   1,   0,   4,     -0.289,     -1.59,   -0.0053 },
  494.    { 0,   1,   0,   2,    -24.420,    -25.10,   -0.3000 },
  495.    { 0,   1,   0,   0,   -668.146,   -126.98,   -0.3997 },
  496.    { 0,   1,   0,  -2,   -165.145,   -165.06,    1.9178 },
  497.    { 0,   1,   0,  -4,     -1.877,     -6.46,    0.0339 },
  498.    { 0,   0,   0,   3,      0.403,     -4.01,    0.0023 },
  499.    { 0,   0,   0,   1,   -125.154,   -112.79,   -0.9781 },
  500.    { 2,   0,   0,   4,      0.213,      1.02,    0.0054 },
  501.    { 2,   0,   0,   2,     14.387,     14.78,    0.2833 },
  502.    { 2,   0,   0,   0,    769.016,    767.96,   10.1657 },
  503.    { 2,   0,   0,  -2,   -211.656,   -152.53,   -0.3039 },
  504.    { 2,   0,   0,  -4,    -30.773,    -34.07,    0.3722 },
  505.    { 2,   0,   0,  -6,     -0.570,     -1.40,    0.0109 },
  506.    { 1,   1,   0,   2,     -2.921,    -11.75,   -0.0484 },
  507.    { 1,   1,   0,   0,   -109.673,   -115.18,   -0.9490 },
  508.    { 1,   1,   0,  -2,   -205.962,   -182.36,    1.4437 },
  509.    { 1,   1,   0,  -4,     -4.391,     -9.66,    0.0673 },
  510.    { 1,  -1,   0,   4,      0.283,      1.53,    0.0060 },
  511.    { 1,  -1,   0,   2,     14.577,     31.70,    0.2302 },
  512.    { 1,  -1,   0,   0,    147.687,    138.76,    1.1528 },
  513.    { 1,  -1,   0,  -2,     28.475,     23.59,   -0.2257 },
  514.    { 1,  -1,   0,  -4,      0.636,      2.27,   -0.0102 },
  515.    { 0,   2,   0,   2,     -0.189,     -1.68,   -0.0028 },
  516.    { 0,   2,   0,   0,     -7.486,     -0.66,   -0.0086 },
  517.    { 0,   2,   0,  -2,     -8.096,    -16.35,    0.0918 },
  518.    { 0,   0,   2,   2,     -5.741,     -0.04,   -0.0009 },
  519.    { 0,   0,   2,   0,   -411.608,     -0.2,    -0.0124 },
  520.    { 0,   0,   2,  -2,    -55.173,    -52.14,   -0.1052 },
  521.    { 0,   0,   2,  -4,      0.025,     -1.67,    0.0031 },
  522.    { 1,   0,   0,   1,     -8.466,    -13.51,   -0.1093 },
  523.    { 1,   0,   0,  -1,     18.609,      3.59,    0.0118 },
  524.    { 1,   0,   0,  -3,      3.215,      5.44,   -0.0386 },
  525.    { 0,   1,   0,   1,     18.023,     17.93,    0.1494 },
  526.    { 0,   1,   0,  -1,      0.560,      0.32,   -0.0037 },
  527.    { 3,   0,   0,   2,      1.060,      2.96,    0.0243 },
  528.    { 3,   0,   0,   0,     36.124,     50.64,    0.6215 },
  529.    { 3,   0,   0,  -2,    -13.193,    -16.40,   -0.1187 },
  530.    { 3,   0,   0,  -4,     -1.187,     -0.74,    0.0074 },
  531.    { 3,   0,   0,  -6,     -0.293,     -0.31,    0.0046 },
  532.    { 2,   1,   0,   2,     -0.290,     -1.45,   -0.0051 },
  533.    { 2,   1,   0,   0,     -7.649,    -10.56,   -0.1038 },
  534.    { 2,   1,   0,  -2,     -8.627,     -7.59,   -0.0192 },
  535.    { 2,   1,   0,  -4,     -2.740,     -2.54,    0.0324 },
  536.    { 2,  -1,   0,   2,      1.181,      3.32,    0.0213 },
  537.    { 2,  -1,   0,   0,      9.703,     11.67,    0.1268 },
  538.    { 2,  -1,   0,  -2,     -2.494,     -1.17,   -0.0017 },
  539.    { 2,  -1,   0,  -4,      0.360,      0.20,   -0.0043 },
  540.    { 1,   2,   0,   0,     -1.167,     -1.25,   -0.0106 },
  541.    { 1,   2,   0,  -2,     -7.412,     -6.12,    0.0484 },
  542.    { 1,   2,   0,  -4,     -0.311,     -0.65,    0.0044 },
  543.    { 1,  -2,   0,   2,      0.757,      1.82,    0.0112 },
  544.    { 1,  -2,   0,   0,      2.580,      2.32,    0.0196 },
  545.    { 1,  -2,   0,  -2,      2.533,      2.40,   -0.0212 },
  546.    { 0,   3,   0,  -2,     -0.344,     -0.57,    0.0036 },
  547.    { 1,   0,   2,   2,     -0.992,     -0.02,    0      },
  548.    { 1,   0,   2,   0,    -45.099,     -0.02,   -0.0010 },
  549.    { 1,   0,   2,  -2,     -0.179,     -9.52,   -0.0833 },
  550.    { 1,   0,  -2,   2,     -6.382,     -3.37,   -0.0481 },
  551.    { 1,   0,  -2,   0,     39.528,     85.13,   -0.7136 },
  552.    { 1,   0,  -2,  -2,      9.366,      0.71,   -0.0112 },
  553.    { 0,   1,   2,   0,      0.415,      0.10,    0.0013 },
  554.    { 0,   1,   2,  -2,     -2.152,     -2.26,   -0.0066 },
  555.    { 0,   1,  -2,   2,     -1.440,     -1.30,    0.0014 },
  556.    { 0,   1,  -2,  -2,      0.384,      0.0,     0.0    },
  557.    { 2,   0,   0,   1,     -0.586,     -1.20,   -0.0100 },
  558.    { 2,   0,   0,  -1,      1.750,      2.01,    0.0155 },
  559.    { 2,   0,   0,  -3,      1.225,      0.91,   -0.0088 },
  560.    { 1,   1,   0,   1,      1.267,      1.52,    0.0164 },
  561.    { 1,  -1,   0,  -1,     -1.089,      0.55,    0      },
  562.    { 0,   0,   2,  -1,      0.584,      8.84,    0.0071 },
  563.    { 4,   0,   0,   0,      1.938,      3.60,    0.0401 },
  564.    { 4,   0,   0,  -2,     -0.952,     -1.58,   -0.0130 },
  565.    { 3,   1,   0,   0,     -0.551,      0.94,   -0.0097 },
  566.    { 3,   1,   0,  -2,     -0.482,     -0.57,   -0.0045 },
  567.    { 3,  -1,   0,   0,      0.681,      0.96,    0.0115 },
  568.    { 2,   0,   2,   0,     -3.996,      0,       0.0004 },
  569.    { 2,   0,   2,  -2,      0.557,     -0.75,   -0.0090 },
  570.    { 2,   0,  -2,   2,     -0.459,     -0.38,   -0.0053 },
  571.    { 2,   0,  -2,   0,     -1.298,      0.74,    0.0004 },
  572.    { 2,   0,  -2,  -2,      0.538,      1.14,   -0.0141 },
  573.    { 1,   1,  -2,  -2,      0.426,      0.07,   -0.0006 },
  574.    { 1,  -1,   2,   0,     -0.304,      0.03,    0.0003 },
  575.    { 1,  -1,  -2,   2,     -0.372,     -0.19,   -0.0027 },
  576.    { 0,   0,   4,   0,      0.418,      0,       0      },
  577.    { 2,  -1,   0,  -1,     -0.352,     -0.37,   -0.0028 }
  578. };
  579.  
  580. # if MOON_TEST_CORR
  581.  
  582. /* moon additional correction terms */
  583.  
  584. struct m5dat {
  585.   REAL8 lng;
  586.   int  i0,i1,i2,i3;
  587. } m5[] = {
  588.  
  589.     /*  lng     l   l'   F    D   */
  590.  
  591.      {  0.127,  0,  0,   0,   6 },
  592.      { -0.151,  0,  2,   0,  -4 },
  593.      { -0.085,  0,  0,   2,   4 },
  594.      {  0.150,  0,  1,   0,   3 },
  595.      { -0.091,  2,  1,   0,  -6 },
  596.      { -0.103,  0,  3,   0,   0 },
  597.      { -0.301,  1,  0,   2,  -4 },
  598.      {  0.202,  1,  0,  -2,  -4 },
  599.      {  0.137,  1,  1,   0,  -1 },
  600.      {  0.233,  1,  1,   0,  -3 },
  601.      { -0.122,  1, -1,   0,   1 },
  602.      { -0.276,  1, -1,   0,  -3 },
  603.      {  0.255,  0,  0,   2,   1 },
  604.      {  0.254,  0,  0,   2,  -3 },
  605.      { -0.100,  3,  1,   0,  -4 },
  606.      { -0.183,  3, -1,   0,  -2 },
  607.      { -0.297,  2,  2,   0,  -2 },
  608.      { -0.161,  2,  2,   0,  -4 },
  609.      {  0.197,  2, -2,   0,   0 },
  610.      {  0.254,  2, -2,   0,  -2 },
  611.      { -0.250,  1,  3,   0,  -2 },
  612.      { -0.123,  2,  0,   2,   2 },
  613.      {  0.173,  2,  0,  -2,  -4 },
  614.      {  0.263,  1,  1,   2,   0 },
  615.      {  0.130,  3,  0,   0,  -1 },
  616.      {  0.113,  5,  0,   0,   0 },
  617.      {  0.092,  3,  0,   2,  -2 },
  618.      {  0.000, 99,  0,   0,   0 }    /* end mark */
  619. };
  620. # endif /* MOON_TEST_CORR */
  621.  
  622.  
  623. /*****************************************************
  624. $Header: deltat.c,v 1.10 93/01/27 14:37:06 alois Exp $
  625. deltat.c
  626. deltat(t): returns delta t (in julian days) from universal time t
  627. is included by users
  628. ET = UT +  deltat
  629.  
  630. ---------------------------------------------------------------
  631. | Copyright Astrodienst Zurich AG and Alois Treindl, 1989.    |
  632. | The use of this source code is subject to regulations made  |
  633. | by Astrodienst Zurich. The code is NOT in the public domain.|
  634. |                                                             |
  635. | This copyright notice must not be changed or removed        |
  636. | by any user of this program.                                |
  637. ---------------------------------------------------------------
  638.  
  639. ******************************************************/
  640.  
  641. double deltat(jd_ad)
  642. double jd_ad;
  643. { 
  644.   static short int ARR dt[] = { /* in centiseconds */
  645.  
  646.     /*
  647.     ** dt from 1637 to 2000, as tabulated in A.E.
  648.     ** the values 1620 - 1636 are not taken, as they fit
  649.     ** badly the parabola 25.5 t*t for the next range. The
  650.     ** best crossing point to switch to the parabola is
  651.     ** 1637, where we have fitted the value for continuity
  652.     */
  653.                 6780, 6500, 6300,
  654.     6200, 6000, 5800, 5700, 5500,
  655.     5400, 5300, 5100, 5000, 4900,
  656.     4800, 4700, 4600, 4500, 4400,
  657.     4300, 4200, 4100, 4000, 3800, /* 1655 - 59 */
  658.     3700, 3600, 3500, 3400, 3300,
  659.     3200, 3100, 3000, 2800, 2700,
  660.     2600, 2500, 2400, 2300, 2200,
  661.     2100, 2000, 1900, 1800, 1700,
  662.     1600, 1500, 1400, 1400, 1300,
  663.     1200, 1200, 1100, 1100, 1000,
  664.     1000, 1000,  900,  900,  900,
  665.      900,  900,  900,  900,  900,
  666.      900,  900,  900,  900,  900,  /* 1700 - 1704 */
  667.      900,  900,  900, 1000, 1000,
  668.     1000, 1000, 1000, 1000, 1000,
  669.     1000, 1000, 1100, 1100, 1100,
  670.     1100, 1100, 1100, 1100, 1100,
  671.     1100, 1100, 1100, 1100, 1100,
  672.     1100, 1100, 1100, 1100, 1200, /* 1730 - 1734 */
  673.     1200, 1200, 1200, 1200, 1200,
  674.     1200, 1200, 1200, 1200, 1300,
  675.     1300, 1300, 1300, 1300, 1300,
  676.     1300, 1400, 1400, 1400, 1400,
  677.     1400, 1400, 1400, 1500, 1500,
  678.     1500, 1500, 1500, 1500, 1500, /* 1760 - 1764 */
  679.     1600, 1600, 1600, 1600, 1600,
  680.     1600, 1600, 1600, 1600, 1600,
  681.     1700, 1700, 1700, 1700, 1700,
  682.     1700, 1700, 1700, 1700, 1700,
  683.     1700, 1700, 1700, 1700, 1700,
  684.     1700, 1700, 1600, 1600, 1600, /* 1790 - 1794 */
  685.     1600, 1500, 1500, 1400, 1400,
  686.     1370, 1340, 1310, 1290, 1270, /* 1800 - 1804 */
  687.     1260, 1250, 1250, 1250, 1250,
  688.     1250, 1250, 1250, 1250, 1250,
  689.     1250, 1250, 1240, 1230, 1220,
  690.     1200, 1170, 1140, 1110, 1060,
  691.     1020,  960,  910,  860,  800,
  692.      750,  700,  660,  630,  600,  /* 1830 - 1834 */
  693.      580,  570,  560,  560,  560,
  694.      570,  580,  590,  610,  620,
  695.      630,  650,  660,  680,  690,
  696.      710,  720,  730,  740,  750,
  697.      760,  770,  770,  780,  780,
  698.      788,  782,  754,  697,  640,  /* 1860 - 1864 */
  699.      602,  541,  410,  292,  182,
  700.      161,   10, -102, -128, -269,
  701.     -324, -364, -454, -471, -511,
  702.     -540, -542, -520, -546, -546,
  703.     -579, -563, -564, -580, -566,
  704.     -587, -601, -619, -664, -644, /* 1890 - 1894 */
  705.     -647, -609, -576, -466, -374,
  706.     -272, -154,   -2,  124,  264,
  707.      386,  537,  614,  775,  913,
  708.     1046, 1153, 1336, 1465, 1601,
  709.     1720, 1824, 1906, 2025, 2095,
  710.     2116, 2225, 2241, 2303, 2349, /* 1920 - 1924 */
  711.     2362, 2386, 2449, 2434, 2408,
  712.     2402, 2400, 2387, 2395, 2386,
  713.     2393, 2373, 2392, 2396, 2402,
  714.     2433, 2483, 2530, 2570, 2624,
  715.     2677, 2728, 2778, 2825, 2871,
  716.     2915, 2957, 2997, 3036, 3072, /* 1950 - 1954 */
  717.     3107, 3135, 3168, 3218, 3268,
  718.     3315, 3359, 3400, 3447, 3503,
  719.     3573, 3654, 3743, 3829, 3920,
  720.     4018, 4117, 4223, 4337, 4449,
  721.     4548, 4646, 4752, 4853, 4959,
  722.     5054, 5138, 5217, 5296, 5379, /* 1980 - 1984 */
  723.     5434, 5487, 5532, 5582, 5630, /* 1985 - 89 from AE 1991 */
  724.     5686, 5757, 5900, 5900, 6000, /* AE 1993 and extrapol */
  725.     6050, 6100, 6150, 6200, 6250, /* 1995 - 1999 */
  726.     6300                          /* 2000 */
  727.   };
  728.  
  729.   double yr, cy, delta;
  730.   long iyr, i;
  731.   yr = (jd_ad + 18262) / 365.25 + 100.0;    /*  year  relative 1800 */
  732.   cy = yr / 100;
  733.   iyr =  (long) (RFloor(yr) + 1800);   /* truncated to integer, rel 0 */
  734.  
  735. #if TIDAL_26    /* Stephenson formula only when 26" tidal term in lunar motion */
  736.  
  737.   if (iyr >= 1637  && iyr < 2000)
  738.   {
  739.     i = iyr - 1637;
  740.     delta = dt[i] * 0.01 + (dt[i+1] - dt[i]) * (yr - RFloor(yr)) * 0.01;
  741.   }
  742.  
  743.   else if (iyr >= 2000)  /* parabola, fitted at value[2000] */
  744.   {
  745.     delta = 25.5 * cy * cy  - 25.5 * 4 + 63.00;
  746.   }
  747.  
  748.   else if (iyr >= 948)   /* from 948 - 1637 use parabola */
  749.   {
  750.     delta = 25.5 * cy * cy;
  751.   }
  752.  
  753.   else   /* before 984 use other parabola */
  754.   {
  755.     delta = 1361.7  + 320 * cy + 44.3 * cy * cy;  /* fits at 948 */
  756.   }
  757.  
  758. #else    /* use Clemence value + 5 sec before 1690, new dt afterwards */
  759.  
  760.   cy -= 1;  /* epoch 1900 */
  761.   if (iyr >= 1690  && iyr < 2000)
  762.   {
  763.     i = iyr - 1637;
  764.     delta = dt[i] * 0.01 + (dt[i+1] - dt[i]) * (yr - RFloor(yr)) * 0.01;
  765.   } 
  766.  
  767.   else if (iyr >= 2000)  /* parabola, fitted at value[2000] */
  768.   {
  769.     delta = 29.949 * cy * cy  - 29.949 * 4 + 63.0;
  770.   } 
  771.  
  772.   else
  773.   {
  774.     delta = 5 + 24.349 + 72.3165 * cy + 29.949 * cy * cy; /* fits at 1690 */
  775.   }
  776.  
  777. #endif /* TIDAL_26 */
  778.  
  779.   return delta / 86400.0;
  780. }
  781.  
  782.  
  783. /*******************************************
  784. $Header: d2l.c,v 1.9 91/11/16 16:24:20 alois Exp $
  785.  
  786. double to long with rounding, no overflow check
  787. *************************************/ 
  788.  
  789. long d2l(x)
  790. double x;
  791. {
  792.   if (x >=0)
  793.     return ((long) (x + 0.5));
  794.  
  795.   else
  796.     return (-(long) (0.5 - x));
  797. }
  798.  
  799.  
  800. /*
  801.  * $Header$
  802.  *
  803.  * A collection of useful functions for centisec calculations.
  804.  
  805. ---------------------------------------------------------------
  806. | Copyright Astrodienst Zurich AG and Alois Treindl, 1991.    |
  807. | The use of this source code is subject to regulations made  |
  808. | by Astrodienst Zurich. The code is NOT in the public domain.|
  809. |                                                             |
  810. | This copyright notice must not be changed or removed        |
  811. | by any user of this program.                                |
  812. ---------------------------------------------------------------
  813. *******************************************************/
  814.  
  815. double degnorm(p)
  816. double p;
  817. {
  818.   if (p < 0) do {
  819.     p += 360.0;
  820.   } while (p < 0);
  821.  
  822.   else if (p >= 360.0) do {
  823.     p -= 360.0;
  824.   } while (p >= 360.0);
  825.  
  826.   return (p);
  827. }
  828.  
  829.  
  830. /*********************************************************
  831. $Header: julday.c,v 1.9 91/11/16 16:25:06 alois Exp $
  832. *********************************************************/
  833.  
  834. /*
  835. ** This function returns the absolute Julian day number (JD)
  836. ** for a given calendar date.
  837. ** The aruguments are a calendar date: day, month, year as integers,
  838. ** hour as double with decimal fraction.
  839. ** If gregflag = 1, Gregorian calendar is assumed, gregflag = 0
  840. ** Julian calendar is assumed.
  841. **
  842. ** The Julian day number is system of numbering all days continously
  843. ** within the time range of known human history. It should be familiar
  844. ** for every astrological or astronomical programmer. The time variable
  845. ** in astronomical theories is usually expressed in Julian days or
  846. ** Julian centuries (36525 days per century) relative to some start day;
  847. ** the start day is called 'the epoch'.
  848. ** The Julian day number is a double representing the number of
  849. ** days since JD = 0.0 on 1 Jan -4712, 12:00 noon.
  850. ** Midnight has always a JD with fraction .5, because traditionally
  851. ** the astronomical day started at noon.
  852. **
  853. ** NOTE: The Julian day number is named after the monk Julianus. It must
  854. ** not be confused with the Julian calendar system, which is named after
  855. ** Julius Cesar, the Roman politician who introduced this calendar.
  856. **
  857. ** Original author: Marc Pottenger, Los Angeles.
  858. ** with bug fix for year < -4711   15-aug-88 by Alois Treindl
  859. **
  860. ** References: Oliver Montenbruck, Grundlagen der Ephemeridenrechnung,
  861. **             Verlag Sterne und Weltraum (1987), p.49 ff
  862. **
  863. ** related functions: revjul() reverse Julian day number: compute the
  864. **             calendar date from a given JD
  865. */
  866.  
  867. double julday(month, day, year, hour, gregflag)
  868. int month;
  869. int day;
  870. int year;
  871. double hour;
  872. int gregflag;
  873. {
  874.   double jd, u, u0, u1, u2;
  875.  
  876.   u = year;
  877.  
  878.   if (month < 3)
  879.     u -=1;
  880.  
  881.   u0 = u + 4712.0;
  882.   u1 = month + 1.0;
  883.  
  884.   if (u1 < 4)
  885.     u1 += 12.0;
  886.  
  887.   jd = RFloor(u0*365.25)
  888.      + RFloor(30.6*u1+0.000001)
  889.      + day + hour/24.0 - 63.5;
  890.  
  891.   if (gregflag) 
  892.   {
  893.     u2 = RFloor(ABS8(u) / 100) - RFloor(ABS8(u) / 400);
  894.  
  895.     if (u < 0.0)
  896.       u2 = -u2;
  897.  
  898.     jd = jd - u2 + 2;
  899.  
  900.     if ((u < 0.0) && (u/100 == RFloor(u/100)) && (u/400 != RFloor(u/400)))
  901.       jd -= 1;
  902.   }
  903.   return jd;
  904. }
  905.  
  906.  
  907. /*********************************************************
  908. $Header: revjul.c,v 1.9 91/11/16 16:25:37 alois Exp $
  909. *********************************************************/
  910.  
  911. /*
  912. ** revjul() is the inverse function to julday(), see the description there.
  913. ** Arguments are julian day number, calendar flag (0=julian, 1=gregorian)
  914. ** return values are the calendar day, month, year and the hour of
  915. ** the day with decimal fraction (0 .. 23.999999).
  916. **
  917. ** Original author Mark Pottenger, Los Angeles.
  918. ** with bug fix for year < -4711 16-aug-88 Alois Treindl
  919. */
  920.  
  921. void revjul(jd, gregflag, jmon, jday, jyear, jut)
  922. double jd;
  923. int gregflag;
  924. int *jmon;
  925. int *jday;
  926. int *jyear;
  927. double *jut;
  928. {
  929.   double u0, u1, u2, u3, u4;
  930.  
  931.   u0 = jd + 32082.5;
  932.   if (gregflag)
  933.   {
  934.     u1 = u0 + RFloor(u0/36525.0) - RFloor(u0/146100.0) - 38.0;
  935.  
  936.     if (jd >= 1830691.5) u1 +=1;
  937.       u0 = u0 + RFloor(u1/36525.0) - RFloor(u1/146100.0) - 38.0;
  938.   }
  939.  
  940.   u2 = RFloor(u0 + 123.0);
  941.   u3 = RFloor((u2 - 122.2) / 365.25);
  942.   u4 = RFloor((u2 - RFloor(365.25 * u3)) / 30.6001);
  943.  
  944.   *jmon = (int)(u4-1.0);
  945.  
  946.   if (*jmon > 12)
  947.     *jmon -= 12;
  948.  
  949.   *jday = (int)(u2 - RFloor(365.25 * u3) - RFloor(30.6001 * u4));
  950.   *jyear = (int)(u3 + RFloor((u4 - 2.0) / 12.0) - 4800);
  951.   *jut = (jd - RFloor(jd + 0.5) + 0.5) * 24.0;
  952. }
  953. #endif /* PLACALC */
  954.  
  955. /* placalc2.c */
  956.